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ABSTRACT 



ey. Lahtornia 93940 



Handling nonlinear systems is not a simple task, 
especially if a computer is not available. Generally, 
although a precise result for every case would be preferable, 
an approximate solution is usually sufficient for most 
engineering solutions. 

R. Boxer and S. Thaler developed a technique [1] 
presented in Proceedings of IRE that allow to obtain the 

K 

response of linear and nonlinear systems, without knowledge 
of the roots of the system characteristic equation. 

The aim of this work is to test this technique in 
obtaining the solution of several time varying linear 
equations and to produce a computer program that improves 
the solution accuracy and ease of use. 
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I . INTRODUCTION 



Two basic methods exist for solving the response of a 
system which generally Is presented as an Integro-dlf ferentlal 
equation: analytic methods and numerical methods. 

For simple problems, analytic methods give precise 
solutions and have the advantage that these solutions can 
be given as a function of the parameters of the system. 

As the order of the differential equation Increases, 
numerical methods give less labor In providing a solution 
and can be, at some extent the only method available. 

Laplace transform Is limited to the case of linear 
differential equations. Nonlinear differential equations 
are out of the scope of this technique, preventing Its 
application. 

EVen sometimes with linear systems. If the Inverse 
Integral Is complex, it Is difficult and sometimes Impossible 
to find analytically the time function. 

However, It Is possible to approximate the Laplace 
transform of a linear system by a z-transform, obtaining 
thereafter a time series representing the values of the 
response at equally spaced intervals of time, by synthetic 
division . 

The method [1] described and used in this paper utilizes 
the z-transform directly in order to obtain the solution of 
the linear and nonlinear system.s. 
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In this paper, an explanation is given of this method 
applied to linear systems and its extension to time-varying 
and nonlinear systems. 

A computer was used for faster and more reliable results, 
but the method described is compatible with handwork. 

C.P.U. time in the computer, for this technique, is 
much smaller than other numerical methods due to the less 
number of iteration steps needed. This technique could 
reduce solution times by several orders of magnitude over 
conventional numerical methods. 



7 



II. PRESENTATION OF THE METHOD 



A. PROCEDURE 

The basic procedures of the method are given, before 
the description or the theory, to show its simplicity. 

Suppose that one has the Laplace transform of an output 
response. It will be usually in the following form: 

a + a, s + a^s^ + • • • + a, s^ + . . . + a„s^ 

PCs) = ^ i ^ ^ (1) 

b +bTS + b„s‘^+ • . . + b, s^ + ... + b s 

o 1 2 k m 



The steps to obtain the time response are: 

1. Express the function P(s) as a rational fraction in 
powers of s“^ by dividing the numerator and denominator 
by s"^. 

-k - 1 

2. Substitute for s a rational fraction in powers of z 

obtained from Table I and rearrange F(s) as a rational 
fraction in powers of z~^. 

3. Divide the resulting expression by T, where T is the 
time Interval between points at which the result is desired. 

4. Expand the fraction by synthetic division into a series 
of the form: 

Dq + D^ z ^ + D 2 z ^ + • • • + D^ z ^ + • • • ( 2 ) 

where D^, the coefficient of z~^, is the approximate value 
of the time response at t = nT. 
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B. LINEAR SYSTEMS 



Take an arbitrary function f(t) as In Fig. 1(a). Its 
Laplace transform and Inverse Laplace transform are defined 
as : 



P(s) = / f(t) e“®^ dt 



(3) 



f(t) = 



2itJ 



C+J" 

/ f(s) e^“ ds 



c-j« 



(4) 





Fig. 1(b) : Sampled f(t) 



Figure 1(b) shows f(t) sampled by shifted unit Impulses. 
The result Is a series of Impulses at t “ nT. 
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Applying Laplace to the sampled function one gets; 



F*(e®^) = E f(nt) (5) 

n=0 

s T 

substituting In (5) e by z so; 

z = e (6) 

the result obtained Is known as the z-transform; 

oo 

F*(z) = Z f(nT) z~^ (7) 

n=0 



and the Inverse Is given by; 

f'(nT) = ifi F*(z) z^”^ dz 



( 8 ) 



where the contour Is a circle centered at the origin and 
Including all the singularities of F*(z). 

The equation (6) maps the jw axis of the s-plane into a 
unit circle around the origin In the z-plane. 

There are other methods for obtaining direct and inverse 
z-transforms . 

In closed form the function will appear as: 



F*(z) 



-1 -2 

a^ + a, z + a^z + ••• + a„z 
o 1 2 ^ 

b« + b, z~^ + b-z”^ + ... + b_z 
0 12 m 



( 9 ) 
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Division of the numerator by the denominator, will 
result in the coefficients of f(nT) as defined in ( 7 ), 
without the use of the inversion integral formula (8). 

In order to apply this useful characteristic of 
z-transform theory to non-sampled systems one uses the 
inverse Laplace integral (^), 

Q+ J C30 

^ c-j“ 



which can be expressed piecewise as : 



f(t) 



, Jtt/T J“ 

^ / f(s)e®'^ds + ^ / [P(s)e®^ + f(-s)e"®^]ds 

-Jtt/T Jtt/T 

( 10 ) 



The way the contour is divided is shown in Pig. 2. The 
0 

first integral gives the part of the contour shown between 
“(tt/T) £ w _< tt/T , The remaining contour u > (tt/T),w < -(tt/T) 
is given by the second integral. 

If 1/T is sufficiently large, a good approximation will 
be held by the first Integral of (10) and the second 
Integral may be discarded, being the resulting error given 
by: 



1 

2ttJ 



J« 

/ 

-Jtt/T 



f(s)e®^ + P(-s)e 



e(t) 



ds 



( 11 ) 




Figure 2; Contour in s-plane for approximate solution 



The approximation, f. (t), of f(t) is given by 






Jtt/T 



St 



/ P(s) ds 



-Jtt/T 



( 12 ) 



Making t = nT 






snT 

/ P(s) ds 

-Ju/T 



(13) 



Solving eq . (6) for s, one obtains: 



s = 1/T In z 



(liJ) 
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Making the substitution of (1^1) into (13) the following 
equation is obtained 



f^(nT) 



1 Jtt/T 
2TTj / 

-Jtt/T 



F(l/T In z) ^ d(l/T In z) 

(15) 



f^(nT) = 2^ ^ P(l/T In z) z^ ^ dz 



(16) 



where the contour, again, is the unit circle in the z-plane. 

If eq . (l6) is compared with eq. (8), it is seen that it 
has the form of an inverse z-transform except for the constant 
multiplier 1/T. 

On the other hand, as F(l/T In z) is a transcendental 
function, it is not expandable directly into a series of 
powers of by synthetic division. The synthetic division 
technique can, however, be used if a proper approximation 
could be found. 

If F(s) is expanded in a series of descending powers of 
s like 

F(s) = a^s“^ + + ••• (17) 

a solution of the form 

^2 

f(t) = a^ + a2t + a^ ^ + ••• (l8) 



can be obtained, but this solution is not always valid. 
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The way to follow is to rearrange P (l/T)/ln z as a 
rational function of two polynomials in powers of s~^ = T/ln 
For this purpose, the power series for In z used is: 

In z = 2[u + 1/3-u^ + 1/5-u^ + •••] (19) 



where : 



u = 



1 - z 
1 + z 



-1 



So the following equation is obtained: 



1 ^ T ^ 

^ ^ ' u + 1/3 u3 + 1/5 u5 



( 20 ) 



By synthetic division the following Laurent series is 
obtained : 



I = T = T r 1 1 „ 3 „5 . 

s Inz 2^n~3 '^ 'PS 



( 21 ) 



To obtain series expansion of s bothsides of (21) 

t h 

can be raised to the power. 

Table I (on page 15) was obtained by using only the 
principal part and the constant term of the Laurent series. 



The result is of the form 



TABLE I 



s-k 


z-Form = Fj^(z 


-1 


T 1 + z"^ 


s 


2 T -1 




1 - z 


-2 


2 -1 -2 
1 + lOz -^ + z 


s 


(l-.-l)2 




z"^ + z“^ 


s 


2 (l-z-l)3 




z"^ + 4z”^ + z“^ T^ 


s » 


T 


c-5 


z“^ + llz”^ + llz”^ + z”^ 


s 


51T (l-z-l)5 
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( 22 ) 



-k NIC . (^-1) 

(1-z ^ 



where Nj^ is a polynomial in powers of z The right hand 
side of (24) is called the "z-form" for s”^. 



Example 1 — Third Order System (Pig. 3) 









0.-- 






Fig. 3: Typical Feedback System for Example 1 



With a step input, as shov/n, the Laplace transform of 
the output is: 



6o(s) = 3 



2 

+ s + s 



(23) 



Following the step by step procedure as given in the 
Introduction , 



6o(s) = 



-3 



1 + s“^ + s“^ 



(24) 



l6 



— k 

Now, substituting the correspondent forms of s from 
Table I and dividing the result by T; 






6T^(z ^ + z~^) 

(12+6Ttr^) - (36+6T-9 T^)z”^ + (36-6T-9T^)z”^ - (12-61H-T^)z“^ 



(25) 



The solution is obtained by choosing T, the time interval 
at which rate our solution is obtained, and proceeding in 
dividing the denominator into the numerator. 

Let us assume T = 1.0, resulting in the following output: 






6z”^ + 6z“^ 



19 - 33z~^ + 21z"^ - 7z"^ 



( 26 ) 



By long division, the result is: 



. 316z ^ + . 86i|z ^ + 1.15z ^ + I.l6z ^ + 1.06z ^ + .987z "^ + • • • 

The accuracy improves if smaller value for T is chosen. 
Letting T = 0.5 in this same example, substituting on eq. (25) 



0 *(z) = l‘5z 

15.25 - 36.75z"^ + 30.75z"^ - 9-25z"^ 



(27) 



carrying out the long division, the result is: 



0 . 984z“^ + . 335 z“^ + . 6 l 0 z“^ + . 853z"^ + • • • 



The points obtained are plotted together with the precise 
solution in Fig. A. 
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X-scale = 2.0 units/inch 
Y-scale = .5 units/inch 



— Precise solution 
X Approximate solution 
for At = 1.0 

+ Approximate solution 
for At = 0.5 



Example 2 — Solution of Linear Differential Equation 

To use the present technique in obtaining the solution 
of integro-dif ferential equations, one takes the Laplace 
transform of the equation, initial conditions included as 
usual, and then proceeds as before. However, after trans- 
forming the equation, the following rules can be applied: 



18 



1. If the degree of the denominator of the Laplace 
transform is higher than the degree of the numerator by two 
or more, proceed directly. 

2. If the degree of the numerator is one less, equal to, 
or greater than the degree of the denominator, divide until 
the remainder meets the criterion of Rule 1. 

Using the procedure stated in the rules above one 
eliminates the high frequency from P(s), that is Implied if 
the numerator is one degree less than the denominator. 

To Illustrate this, a simple example is shown: 



dt 



+ 



Y = 1 



( 28 ) 



Assuming for initial condition Y(0) = u, the Laplace 
transform Y(s) is: 



Y(s) 



l+s*u_s*u+l 
irmr- 32 + 3 



Using Rule 2 above, one gets: 



( 29 ) 



Y(s) = + -- - -- - ( 30 ) 

s^ + s 

In this procedure, the step function of amplitude u is 
added to the result obtained after operating upon 
(1 - u)/(s^ + s) . 
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Now, letting u = 0.1, and expressing the second term 
of eq. ( 30 ) In powers of s“^. 



Y^^Cs) 



0.9s 



-2 



1 + s 



-1 



( 31 ) 



Substituting the z-forms of Table I (on page 15) and 
dividing by T, 






0.9T(1 + lOz ^ + z ^) 



(12 + 6T) - 2ilz“^ + (12-6 T)z ^ 



(32) 



Letting T = 0.1 and carrying out the long division one 
gets : 



.00714 + .085z ^ + .i 63 z + 



-2 



.233z ^ + 



.296z 



-4 



(33) 



This result Is plotted in Fig. 5 after adding the step 



function of amplitude 0.1. 




Oj'i 910 

Figure 5 

X-scale =2.0 units/inch — Precise solution 

Y-scale = .5 units/inch + Approximate solution 
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C. TIME-VARYING SYSTEMS 

A linear differential equation with time-varying 
coefficients may be expressed as 



d 



n 



dt 



n 



fn(t)Y 



d 



n-1 



dt 



n-1 






+ fj(t)Y = f(t) (34) 



where f^(t), f^_^(t ) , . . . jf^Ct) are functions of the 
independent variable t. 

To solve this type of equation the following steps may 
be executed: 

1. Obtain the Laplace transform of (3^) in the usual 
manner including initial conditions. 

2. Proceed as in the constant coefficient case discussed 
above. In performing the long division process, hov;evsr, 
change the values of C^,...,Cq at each step in the division 
process to the values at the corresponding times of the 
functions 



^n = 



C„ = 



fn(t) 



fo(t) 



(35) 



In general 



C, = 



f^Ct) 



(36) 
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D. NONLINEAR SYSTEMS 



A Similar manner like used in II. C. may be applied to 
nonlinear differential equations. The regression equation 
is obtained by long division process as in the time varying 
example. 

However, the factors , C C^. must now 

represent functions of the dependent variable. 

Generally, a nonlinear differential equation can be 
represented as: 



.n 



dt 



n 



f CY — )Y 
n^ ’dt’ * * * ^ 



.n-1 



dt 



n 



dY 



^n-l^^’dt ’ * • • 



+ = f(t) 



(37) 



being the rules for the procedure as follow; 

1. Obtain the Laplace transform of (37) in the usual 
manner including initial conditions. 

2. Proceed as in the constant coefficient case. In 
performing the long division process, however, change the 
value of C • • • at each step in the division process to 
the most recent values available for 



^n ^n^^’dt’*"^ 



^n-1 ^n-l^^’df * 



( 38 ) 



C = f CY ) 

^o ^o'‘ ’dt’***'^ 
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In general 



Ci = (39) 

In other v;ords, in the first step, the values of are 

H V 

obtained from initial conditions for Y, — , ... . The 
result of the first division is a new value of Y to be used 
in calculating C^, ..., for the second step in the 

division. New values to be used for the derivatives of Y 
with respect to t are calculated from the new values of the 
dependent variables. The process is difficult to describe 
in words, but is simple in concept. 

To illustrate the procedure an example is given: 

Example 3: 



Let 



If + Y^ = 1 

dt 



(^0) 



Carrying out the procedure outlined above, eq.(i|b) is 
written 



If + C-Y = 1 
dt 



(ill) 



The Laplace transform, assuming C to be constant and 
expressing in powers of s”^ is 



Y(s) = ^ 

1 + Cs“-^ 



(i|2) 
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Substituting the z-forms of Table I and dividing by T yields 



Y *(z) = 



T(1 + lOz ^ + z"^) 



(12 + 6CT) - 24z -^ + (12-6CT)z‘ 



Letting T = 0.1, 



Y *(z) = 



0.1 + l.Oz ^ + z ^ 



(12 + 0.60 - 2ilz -^ + (12-0.6C)z‘ 



Carrying out the long division procedure yields 



.00833 + O.lOOz"^ + 0.199z"^ + 0.295z"^ + 0. 390z ^ + 



(^ 3 ) 







I 



SCO IaCC roc 




Figure 6 

X-scale = 2.0 units/inch 
Y-scale = .5 units/inch 



— Precise solution 
+ Approximate solution 
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III. CASE STUDIES 



A. INTRODUCTION 

In order to test the present technique, a collection 
of linear differential equations with varying coefficients 
were studied, whose basic format Is: 



y + P(t)Y = Q(t) 



where P(t) and Q(t) are time functions. For the purpose of 
this work, nine equations In all were processed successfully 
for what P(t) took the values \/t, t and t^, and Q(t) took 
the values l*u(t), t and t . 

The collection of equations, the case nubmers, and the 
sections In this work, are as follows : 



y +ftY = 1 

Y + tY = 1 

Y + t^Y = 1 



Case 1 

Case 2 III B 

Case 3 



Y +/t Y = t 

Y + tY = t 

Y + t^Y = t 



Case 4 

Case 5 III C 

Case 6 



Y +'ftY = 

Y + tY = t^ 

Y + t^Y = t^ 



Case 7 

Case 8 III C 

Case 9 



In Case 10, an example from Ref. 2 v;as studied. Comments 
on this problem are given In Case 10 discussion. 
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In cases 1 to 9 the same type of procedure was followed, 
applying directly the steps as in II A of this work and 
experimenting for time intervals form 0.2 sec. to 2.0 sec. 
in order to show the advantages of this technique with 
relatively large At’s, its stability together with an 
algorithm for computer solution and application to hand 
calculation, l.e., small number of division steps for good 
results. 

B . CASE l-Y+\/tY = 1 

The equation Y + P(t)Y = Q(t) is rewritten for this 
case as : 



dY 



jtt = 



^(t=0) = “ 



following the procedure, as in II A; 



^ + PY = 1 
dt 



P =Jt 



the Laplace transform is 



sy + Py = 1/s 



rearranging and expressing in powers of s 



-1 



s y + Psy = 1 



yCs'^+Ps) = 1 



-2 



y = 



s^ + Ps 1 + Ps ^ 
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— k 

Substituting the z-form values for s from Table I (page 15) 



y^*(z) = 



1 + lOz ^ 

12 (l-z~^)^ 




Dividing by T and rearranging: 



y^*(z) = 



T + 10Tz“^ + Tz 



(^ 5 ) 



(12 + 6PT) - 24z ^ + (12-6PT)z“^ 



where again P is the time varying coefficient of Y in the 
differential equation and T is the iteration time, or in 
other v;ords the time difference between results. 

The following step consists of dividing the denominator 
into the numerator and substituting in P, its value for the 
corresponding time (t = nT) . As the exponent of z for the 
first value in the quotient is zero, P(t) =7t at t = 0 , 

P(t) = 0 in the first division step. 

The value for T first chosen was 0.2 sec. . 

The algorithm used for computing the quotient (Program 1, 
page 30) as the solution, for values of T from 0.2 sec. up 
to 2 sec. with the "precise" solution in continuous line and 
with the approximate solution as crosses follow; the "precise" 
solution was obtained by computer, using the subroutine 
called INTEGl. 
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Scales in the graphs (Figs. 7 to 20 ) are the same in 
each case. 

In the numerical output tables, the time is referred to 
as T, the approximate solution is referred to as Q(N) and 
the "precise" solution as PR. TPIN is the final time. 

It can be seen that the approximated solution even for 
large time iterations is very close to the precise solution. 



FOR DELTA T = 0.20 TFIN = 5.2 



T 


Q(N) 


PR 


0.0 


0.01667 


0.0 


0.2 


0.19144 


0.19300 


0 . it 


0. 36010 


0.36197 


0.6 


0.49870 


0.50039 


0.8 


0.60588 


0.60716 


1.0 


0.68335 


0.68414 


1.2 


0.73455 


0.73486 


1.4 


0.76372 


0.76363 


1.6 


0.77528 


0.77491 


1.8 


0.77345 


0.77289 


2.0 


0.76192 


0.76128 


2.2 


0.74384 


0.74318 


2.4 


0.72170 


0.72109 


2.6 


0.69743 


0.69689 


2.8 


0.67245 


0.67200 


3.0 


0.64773 


0.64738 


3.2 


0.62392 


0.62366 


3.4 


0.60141 


0.60122 


3.6 


0.58039 


0.58027 


3.8 


0.56092 


0.56085 


4.0 


0.54297 


0.54294 


4.2 


0.52648 


0.52647 


4.4 


0.51132 


0.51133 


4.6 


0.49738 


0.49741 


4.8 


0.48454 


0.48458 


5.0 


0.47268 


0. 47272 
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non o o o , o o noon 



cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

Q(T)=STEP , P(T)=SQRT(T) 

ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

INTEGER *4ITB( 12)712*0/ 

REAL *4RTB( 28) 728*0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(100)» PR(IOO) 
ITB(3) = 5 
ITB(4) = 5 

DO 1 1=1,26 
READ (5,5) XX(I) ,PR(I ) 

1 CONTINUE 

2 READ (5,5,END=4) TD,TFIN 
WRITE (6,6) TD,TFIN 
WRITE (6,7) 

M = TFIN/TD 
A(l) = TD 
A(2) = 10.0*TD 
A(3) = TD 
A(4) = 0.0 

DO 3 N=1,M 
T = (N-1)*TD 

P = SQRT(T) 

FI = 12.0+6.0*( P*TD) 

F2 = -24.0 
F3 = 12.0-6.0*(P*TD) 

Q(N) IS THE OUTPUT 

Q(N) = A(N)/F1 
A(N+1) = A(M+1)-(Q(N)*F2) 

A(N+2) = A(N+2)-(Q(N)*F3) 

A(N+3) = 0.0 
X(N) = T 

WRITE (6,5) T,Q(N) 

3 CONTINUE 
C 

WRITE (6,8) 

WRITE (6,5) (XX( I) ,PR(I) ,1=1,26) 

ITB(l) = 1 
ITB(2) = 0 
ITB(12) = 1 

CALL DRAWP ( 26 , XX, PR , I TB, RTB ) 

ITB(l) = 3 
ITB(2) = 2 

CALL DRAWP ( M , X , Q , I TB, RTB ) 

GO TO 2 

4 STOP 



5 

6 

7 

8 



FORMAT 

FORMAT 

FORMAT 

FORMAT 

END 



( 2F10.5) 

( • 1* , • FOR 

(• 



DELTA T=* , F6.2 ,4X, ' TFIN=',F4.1) 



( 



XX = * , 7X, • 



Q(N)=* ) 
PR=* ) 



PROGRAM 1 
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CC2 (X'i ogs OQ'9 




x-scale =1.0 units/inch 
y-scale = 0.2 units/inch 



FOR DELTA T = 0.2 



FIGURE 7 
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OGO l/CO ZZO 




x-scale =1.0 units/inch 



FOR DELTA T = 0.3 



TPIN =5.2 



T 


Q(N) 


0.0 


0.02500 


0.3 


0.27722 


0.6 


0.49673 


0.9 


0.64695 


1.2 


0.73424 


1.5 


0.77180 


1.8 


0.77421 


2.1 


0.75441 


2.4 


0.72252 


2.7 


0.68563 


3.0 


0.64822 


3.3 


0.61282 


3.6 


0.58059 


3.9 


0.55187 


4.2 


0.52653 


t».5 


0.50422 


4.8 


0.48454 



FIGURE 8 
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FOR DELTA T = O.MO TFIN = 5.2 



T 


Q(N) 


0.0 


0.03333 


0.4 


0.35508 


0.8 


0.60241 


1.2 


0.73386 


1.6 


0.77661 


2.0 


0.76404 


2.4 


0.72370 


2.8 


0.67393 


3.2 


0.62484 


3.6 


0.58087 


4.0 


0.54317 


4.4 


0.51137 


4.8 


0.48453 



FIGURE 9 
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FOR DELTA T = 0.5 



T 



0.0 

0.5 

1.0 

1.5 
2.0 

2.5 
3.0 

3.5 
^ 1.0 

4.5 



TFIN =5.2 

Q(N) 

0.4167 

0.42489 

0.67982 

0.77314 

0.76570 

0.71310 

0.64984 

0.59170 

0.54331 

0.50424 



FIGURE 10 
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Ki^Z OQ-; 0Q6 




FOR DELTA T = 



T 



0.0 

0.6 

1.2 

1.8 

2 . 1 ) 

3.0 

3.6 

4.2 



0.6 TFIN = 

Q(N) 

0.05000 
0.48686 
0. 73288 
0.77864 

0.72725 

0.65099 

0.58165 

0.52673 



5 



.2 



FIGURE 11 
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POO 900 . 200 



o 




\Jju 



QQ2 


003 


DELTA T = 0.7 


TFIN = ! 


T 


Q(N) 


0.0 


0.5833 


0.7 


0.54145 


1.4 


0.76577 


2.1 


0.76211 


2.8 


0.67831 


4.2 


0.52680 



FIGURE 12 



36 



UC2 aJ-1 CTJ5 0Q»^ 




FOR DELTA T = 0.8 
T 



0.0 

0.8 

1.6 

2.H 

3.2 

H.O 



TFIN =5.2 

Q(N) 

0.06667 
0.58920 
0 . 78249 
0.73260 
0.62874 
0.54381 



FIGURE 13 
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0C2 OQ'5 OQG OOS 




FOR DELTA T = 0.9 TFIN = 5.2 



T 


Q(N) 


0.0 


0.07500 


0.9 


0.63073 


1. 8 


0.78658 


2.7 


0.69661 


3.6 


0.583^0 



FIGURE 1 ^ 
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UC2 OO; 02S oay 




FOR DELTA T = 1.0 



TFIN +5.2 



T 



Q(N) 



0.0 

1.0 

2.0 

3.0 

4.0 



0.08333 

0.66667 

0.78105 

0.65849 

0.54411 



FIGURE 15 
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SCO SOO 2X 




FOR DELTA T = 1.2 TFIN =5.2 
T Q(N) 



0.0 0.10000 

1.2 0.72408 

2.4 0.75053 

3.6 0.58590 



FIGURE 16 
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PDO 




FOR DELTA T = 1.4 TPIN =5.2 



T 



Q(N) 



0.0 

1.4 

2.8 



0.11667 

0.76576 

0.7053^ 



FIGURE 17 
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’JC2 0Q6 00^3 




FOR DELTA T = 1.6 TPIN = 5.2 



T 


Q(N) 


0.0 


0.13333 


1.6 


0.79526 


3.2 


0.6542^1 



FIGURE 18 
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CCO 900 ’to 20C COO 




FOR DELTA T = 1.8 TPIN =5.2 



T 



0.0 

1.8 



Q(N) 

0.15000 

0.81541 



FIGURE 19 
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'X2 (Xli 0Q6 OOrJ 



» 







r\'y 1 



no7 



OQ'3 



DELTA T = 2.0 TFIN = 5.2 



T 



Q(N) 



0.0 0.16667 

2.0 0.82843 



FIGURE 20 



4 4 



1 . Case 2 - Y + tY = 1 



The calculations are the same as in Case 1 with the 
only difference that P(t) is now equal to t. So, -when 
dividing as in Equation ( 45 ), repeated here for convenience, 
the denominator into the numerator the substitution for each 
step is as P = t . 



Ya'Cz) 



T + lOTz ^ 

(12 + 6 PT) - 24 z“^ 



+ Tz 

+ (12 - 6 PT)z“^ 



In the same manner, numerical and graphical output for 
T from 0.2 sec. to 2.0 sec. (Pigs. 21 to 34 ) follow. 

Program 2 (page 4 ? ) shows a computer algorithm for solution 
of this problem. 



FOR DELTA T = 0.2 TPIN = 5.2 



T 


Q(N) 


PR 


0.0 


0.01667 


0.0 


0.2 


0.19608 


0.19735 


0.4 


0. 37707 


0.37933 


0.6 


0.53018 


0.53293 


0.8 


0.64664 


0.64931 


1.0 


0. 72264 


0.72478 


1.2 


0.75927 


0.76058 


1.4 


0.76154 


0.76198 


1.6 


0.73700 


0.73670 


1.8 


0.69413 


0.69335 


2.0 


0.64099 


0.64000 


2.2 


0.58425 


0.58329 


2.4 


0.52880 


0.52801 


2.6 


0.47769 


0. 47713 


2.8 


0.43241 


0.43208 


3.0 


0.39333 


0.39319 


3.2 


0 . 36010 


0. 36008 


3.4 


0.33199 


0.33204 


3.6 


0.30817 


0.30825 


3.8 


0.28784 


0.28793 
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FOR DELTA T = 0. 


2 TFIN 

(Continued) 


= 5.2 


T 


Q(N) 


PR 


i|.0 


0.27033 


0.27041 


^,2 


0.25507 


0.25514 


H.n 


0.24162 


0.24168 


4.6 


0.22966 


0.22971 


4.8 


0.21893 


0.21897 


5.0 


0.29023 


0.20926 



cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

C C(T)=STEP T P(T)=T 

c 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

INTEGER *4ITB( 12) /12*0/ 

REAL *4RTB(28)/28=i'0.0/ 



DIMENSION A(IOO) 
ITB(3) = 5 
ITB(4) = 5 



Q(100)t X(IOO), XX(100)» PR(IOO) 



DO 1 1=1,26 

READ (5,7) XX(I),PR(I) 

CONTINUE 



READ (5,7,EN0=4) TD,TFIN 
WRITE (6,5) TO,TFIN 
WRITE (6,6) 

M = TFIN/TD 



A( 1) 
A(2) 
A(3) 
A(4) 



TD 
10 . 0 ^ 
TD 
0.0 



TD 



DO 3 N=1,M 
T = (N-1)=?=TD 

P = T 



C 

C 



FI = 12.0+6.0«(P*TD) 

F2 = -24.0 

F3 = 12 .0-6.0*(P#TD) 

Q(N) IS THE OUTPUT 

Q(N) = A(N)/F1 

A(N<-1) = A(N<-1)-(Q(N)=«''F2) 

A(N<-2) = A(N<-2)-( Q(N)*F3) 

A(N+3) = 0.0 
X(N) = T 

WRITE (6,7) T,Q(N) 

3 CONTINUE 

WRITE (6,8) 

WRITE (6,7) (XX( I ) ,PR( I ) ,1=1,26) 
ITB(l) = 1 
ITB(2) = 0 
ITB(12) = 1 

CALL DRAWP ( 26 , XX , PR , I TB, RTB ) 
ITB( 1) = 3 
ITB(2) = 2 

CALL DRAWP ( M , X , Q , I TB , RTB ) 

GO TO 2 

4 STOP 



5 

6 

7 

8 



FORMAT 

FORMAT 

FORMAT 

FORMAT 

END 



DELTA T=' ,F6.2,4X,' 
,7X,' Q(N)=') 

( 2F10.5) 



• 1' , 'FOR 

I T = ' 



TFIN = * ,F4.1) 



( 



XX=' ,7X 



PR=' ) 



PROGRAM 2 
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0C2 OO'i 006 OOd 




x-scale = 1.0 units/inch 
y-scale = 0.2 units/inch 

FOR DELTA T = 0.2 



FIGURE 21 
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UC 2 024 006 COa 




DELTA T = 0.3 


TFIN = ' 


T 


Q(N) 


0.0 


0.02500 


0.3 


0.28708 


0.6 


0.52675 


0.9 


0.68665 


1.2 


0.75759 


1.5 


0.75202 


1.8 


0.69513 


2.1 


0. 61^102 


2.il 


0.52986 


2.7 


0 . 45 'i 88 


3.0 


0.39355 


3.3 


0 . 3 '< 5''5 


3.6 


0.30809 


3.9 


0.27869 


il.2 


0.25500 


' 1.5 


0 . 235''3 


il.8 


0.21891 



FIGURE 22 
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f '20 900 1-00 zon 




DELTA T = 0.40 


TFIN = 


T 


Q(N) 


0.0 


0.03333 


0.4 


0.37037 


0.8 


0.63857 


1.2 


0.75516 


1.6 


0.73782 


2.0 


0.64408 


2.4 


0.53138 


2.8 


0.43354 


3.2 


0. 36022 


3. 6 


0.30795 


4.0 


0.27013 


4.4 


0.24150 


4.8 


0.21887 



FIGURE 23 
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•0t2 OQl 005 00>J 



f » 
O 

o 




FOR DELTA T = 0.5 TFIN =5.2 



T 


Q(N) 


0.0 


0.04167 


0.5 


0.44444 


1.0 


0.71111 


1.5 


0.75151 


2.0 


0.64646 


2.5 


0.50660 


3.0 


0.39427 


3.5 


0.31923 


4.0 


0.26995 


4.5 


0.23529 



FIGURE 24 
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012 OO't 005 ooa 



o 

Cl 

o 




0J2 


0\jZf QO i 


DELTA T = 0.6 


TFIN =5.2 


T 


Q(N) 


0.0 


0.5000 


0.6 


0.50847 


1.2 


0.74776 


1.8 


0.70037 


2.4 


0.53615 


3.0 


0.39480 


3.6 


0.30744 


4.2 


0.25460 



FIGURE 25 
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FOR DELTA T 
T 

0.0 

0.7 

l.ii 

2.1 

2.8 

3.5 

4.2 



= 0.7 



TFIN =5.2 
Q(N) 

0.05833 

0.56225 

0.75470 

0.62530 

0.43722 

0.31854 

0.25438 



FIGURE 26 
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00 1 006 000 




DELTA T = 0.8 


TPIN = ' 


T 


Q(N) 


0.0 


0.06667 


0.8 


0.60606 


1.6 


0.73910 


2.4 


0.54392 


3.2 


0.36042 


4.0 


0.26888 



FIGURE 27 
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fCO 900 TOO 2Xi 




FOR DELTA T = 0.9 
T 



0.0 

0.9 

1.8 

2.7 

3.6 



TFIN =5.2 

• Q(N) 

0.07500 

0.64057 

0.70781 

0.46704 

0.30519 



FIGURE 28 
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PCO 9C0 T'CO 




FOR DELTA T 
T 

0.0 

1.0 

2.0 

3.0 

n.o 



= 1.0 



TFIN = 5.2 

Q(N) 

0.08333 

0.66667 

0.66667 

0. ifOOOO 

0.26667 



FIGURE 29 
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FOR DELTA T = 1.2 



TFIN = 5. 
Q(N) 



T 



0.0 


0.10000 


1.2 


0.69767 


2.H 


0.57186 


3.6 


0.30012 



FIGURE 30 



FOR DELTA T = 1.^1 



TFIN = 5. 
Q(N) 



T 



0.0 


0.11667 


lA 


0.70707 


2.8 


0.^775 



FIGURE 31 
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FOR DELTA T = 1.6 TFIN =5.2 



T 


Q(N) 


0.0 


0.13333 


1.6 


0.70175 


3.2 


0.39^24 



FIGURE 32 
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POO 900 2u0 




FOR DELTA T = 1.8 TFIN = 5.2 



T 



0.0 

1.8 



Q(N) 

0.15000 

0.68702 



FIGURE 33 
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■X2 Qu; OQt; ooa 




FOR DELTA T = 2.0 TFIN =5.2 



T 



Q(N) 



0.0 0.16667 

2.0 0.66667 



FIGURE 34 
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2 . Case 3 - Y + t^Y = 1 

As in Case 1 and Case 2 , the same calculations are 
followed and only in the division, the substitution of 
P(t) = t^ is the sole difference. 

Output for several time iterations follows as in the 
previous cases (Figs. 35 to 48 ). Computer algorithm is shown 
in Program 3 (page 63). 



FOR DELTA T = 0.2 TFIN = 5.2 



T 


Q(N) 




0.0 


0.01667 


0.0 


0.2 


0.19920 


0.19960 


0.4 


0.39213 


0. 39366 


0.6 


0.56550 


0.56858 


0 . 8 


0.70032 


0. 70472 


1.0 


0.77773 


0.78240 


1.2 


0.78667 


0.79011 


1.4 


0.73026 


0.73132 


1.6 


0.62669 


0.62542 


1.8 


0.50321 


0.50080 


2.0 


0.38584 


0.38376 


2.2 


0.29077 


0.28980 


2.4 


0.22210 


0.22208 


2.6 


0.17552 


0.17589 


2.8 


0.14398 


0.14433 


3.0 


0.12163 


0.12183 


3.2 


0.10482 


0.10492 


3.4 


0.09160 


0.09165 


3.6 


0.08088 


0.08091 


3.8 


0.07204 


0.07205 


4.0 


0.06462 


0.06464 


4.2 


0.05833 


0.05834 


4.4 


0.05294 


0.05295 


4.6 


0.04828 


0.04830 


4.8 


0.04422 


0.04423 


5.0 


0.04067 


0.04067 



Equation ( 45 ) is repeated here for convenience: 



y^*(z) 



T + lOTz ^ + Tz~^ 

(12 + 6 PT) - 24 z"^ + (l 2 - 6 PT)z"^ 
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cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

C Q(T)=STEP , P(T)=T**2 
C 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

INTEGER *4ITB( 12)/12*0/ 

REAL ^4RTB(28)/28*0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(IOO), PR(IOO) 
ITB(3) = 5 
ITBI4) = 5 
C 

DO 1 1=1,26 

READ (5,7) XX( I),PR(I ) 

1 CONTINUE 
C 

2 READ (5,7,END=4) TD,TFIN 
WRITE (6,5) TD,TFIN 
WRITE (6,6) 

M = TFIN/TD 
All) = TD 
A(2) = 10.0=f^TD 
A(3) = TD 
A(4) = 0.0 
C 

DO 3 N=1,M 
T = (N-1)*TD 
C 

P = T**2 
C 

FI = 12.0+6.0=4=(P={'TD) 

F2 = -24.0 

F3 = 12.0-6.0«(P*TD) 

C 

C Q(N) IS THE OUTPUT 
C 

Q(N) = A(N)/F1 

A(N+i) = A(N+1)-(Q(N)*F2) 

A(N+2) = A(N+2)-(Q(N)=i=F3) 

A(N+3) = 0.0 
X(N) = T 

WRITE (6,7) T,Q(N) 

3 CONTINUE 
C 

WRITE (6,8) 

WRITE (6,7) (XX(I),PR(I),I=1,26) 

ITB(l) = 1 
ITB(2) = 0 
ITB(12) = 1 

CALL ORAWP ( 26,XX,PR, ITB,RTB) 

ITB(l) = 3 
ITB(2) = 2 

CALL DRAWP ( M , X , Q , I TB , RTB ) 

GO TO 2 

4 STOP 
C 

5 FORMAT C1',*F0R DELTA T= ' , F6 . 2 , 4X , ' TFIN = ',F4.1) 

6 FORMAT (• T=',7X,' Q(N)=') 

7 FORMAT (2F10.5) 

8 FORMAT (• XX=',7X,' PR=') 

END 



PROGRAM 3 
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PCQ &00 too ^cn 




x-scale =1.0 units/inch 
y-scale = 0.2 units/inch 

FOR DELTA T = 0.2 



FIGURE 35 
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PvO 900 1-00 20r. 




FOR DELTA T = 0 . 3 TPIN = 5.2 



T 


Q(N) 


0.0 


0.02500 


0.3 


0.29600 


0.6 


0.56168 


0.9 


0.7^128 


1.2 


0.78225 


1.5 


0.68283 


1.8 


0.50631 


2.1 


0.33719 


2.i< 


0.22218 


2.7 


0.15773 


3.0 


0.12138 


3.3 


0.09778 


3.6 


0.08086 


3.9 


0.06816 


iJ.2 


0.05832 


^.5 


0.05053 


i<.8 


Q . Qi \ i \22 



FIGURE 36 
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UC2 iX>'* 0Q8 00^ 




FOR DELTA T = 0 . TFIN =5.2 



T 


Q(N) 


0.0 


0.03333 


0.4 


0.38760 


0.8 


0.68723 


1.2 


0.77582 


1.6 


0.62989 


2.0 


0.39299 


2.4 


0.22240 


2.8 


0.14260 


3.2 


0.10466 


3.6 


0.08082 


4.0 


0.06460 


4.4 


0.05293 


4.8 


0.04422 



FIGURE 37 
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y02 001 OOS OOiJ 




FOR DELTA T = 0.5 TFIN =5.2 



T 


Q(N) 


0.0 


0.0^1167 


0.5 


0.il7059 


1.0 


0.7529^ 


1.5 


0.68141 


2.0 


0.39906 


2.5 


0.19512 


3.0 


0.12007 


3.5 


0.08613 


H,0 


0.06447 


i*.5 


0.05057 



FIGURE 38 
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FOR DELTA T = 0.6 


TFIN = 


T 


Q(N) 


0.0 


0.05 


0.6 


0.54152 


1.2 


0.75631 


1.8 


0.52210 


2.4 


0.22530 


3.0 


0.11783 


3.6 


0.08177 


4.2 


0.05783 



FIGURE 39 
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0C2 OO-l OQ6 GO'S 




FOR DELTA T = 0.7 TFIN = 5.2 



T 


Q(N) 


0.0 


0.05833 


0.7 


0.59752 


1.4 


0.70881 


2.1 


0.36272 


2.8 


0.13431 


3.5 


0.08809 


4.2 


0:05721 



FIGURE 40 
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0 C 2 


233 C 2 'i 


DELTA T = 0.8 


TFIN = 5.2 


T 


Q(N) 


0.0 


0.06667 


0.8 


0.63694 


1.6 


0.62939 


2.4 


0.23756 


3.2 


0.09620 


4.0 


0.06786 



FIGURE 41 
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UC2 006 QQO 




FOR DELTA T = 0.9 TFIN =5.2 



T 



Q(N) 



0.0 

0.9 

1.8 

2.7 

3.6 



0.07500 

0.65958 

0.53668 

0.15283 

0.08072 



FIGURE 
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0C2 OO'l 000 00 tl 




FOR DELTA T 
T 

0.0 

1.0 

2.0 

3.0 

4.0 



= 1.0 



TPIN =5.2 
Q(N) 

0.08333 

0.66667 

0.44444 

0.10101 

0.07183 



FIGURE 43 
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UC2 00-; ao6 cce 




FOR DELTA T = 1.2 TFIN = 5.2 

T Q(N) 

0.0 0.10000 

1.2 0.64378 

2.4 0.28895 

3.6 0.05587 



FIGURE 44 
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Lvi 032 2C3 03 

FOR DELTA T = l.il TFIN = 5.2 



T 


Q(N) 


0.0 

l.il 

2.8 


0.11667 

0.59022 

0.1819^ 




FIGURE iJ 5 



7^1 



FOR DELTA T = 1.6 



TFIN =5.2 
Q(N) 



T 



0.0 

1.6 

3.2 



0.13333 

0.52493 

0.11422 



FIGURE 46 
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I 

i 



0C2 005 000 OOd 




FOR DELTA T = 1.8 TFIN = 5.2 

T Q(N) 

0.0 0.15000 

1.8 0.^5965 



FIGURE ^7 
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3 

4 



E-CO aDO 




FOR DELTA T = 2.0 TFIN =5.2 



T 



Q(N) 



0.0 0.16667 

2.0 0.40000 



FIGURE 48 



77 



C. CASE 4-Y+VtY = t 

For this case, in the equation Y + P(t)Y = Q(t), 
P(t) = >/t and Q(t) = t, so following the procedure the 
equation is rewritten as; 



^ + PY = t 
at 



^(t=o) ° 



where P = J~t . The Laplace transform is 



sy + Py = 1/s' 



so 



s^y + Ps^y = 1 



y(s^ + Ps^) = 1 



.-3 



y = 



l^+Ps^ 1 



+ Ps 



-1 



Substituting the z-forms of s”^ from Table I one gets 



T^ z”^ + z”^ 



yy^*(z) 



2 d-z"^)^ 



1 + P 



T 1 + z 



-1 



2 1 - z 



-1 



Dividing by T and ordering in crescent powers of z 
in numerator and denominator: 



-1 



T^z ^ + T^z ^ 

y *(z) = ^ ^ 

^ (2+Pt) - (6+PT)z"^ + (6-PT)z“^ + (PT-2)z"3 



(il6) 
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Now, substituting for the time iteration wanted and 
replacing P = >/t for its value at the division step, the 
solution is obtained. 

In the same way as before, plots and outputs follow 
for time iterations from 0.2 sec. to 2.0 sec. (Figs. hS to 
60 ). 

Algorithm for computer solution of this problem is shown 
in Program 4 (page 8o) • 



FOR DELTA T = 0.2 TFIN = 5.2 



T 


Q(N) 


PR 


0.2 


0.01914 


0.0 


0 . 


0.07363 


0.01949 


0.6 


0.15683 


0.07449 


0.8 


0.26131 


0.15817 


1.0 


0.37994 


0.26299 


1 . 2 


0.50647 


0.38181 


l.il 


0.63576 


0.50837 


1.6 


0.76391 


0.63756 


1.8 


0.88813 


0.76551 


2.0 


1.00662 


0.88948 


2.2 


1.11838 


1.00770 


2 .H 


1.22302 


1.11920 


2.6 


1.32060 


1.22360 


2.8 


1.41147 


1.32100 


3.0 


1.49613 


1.41170 


3.2 


1.57520 


1 . 49620 


3 .^ 


1.64928 


1.57520 


3.6 


1.71898 


1.64930 


3.8 


1.78485 


1.71890 


4.0 


1.84739 


1.78480 


4.2 


1.90703 


1.84730 


4.4 


1.96413 


1.90700 


4.6 


2.01902 


1.96410 


4.8 


2.07196 


2.01900 


5.0 


2.12316 


2.07200 


5.2 


2.17281 


2.12320 



79 



f 

\ 



nno 



CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC3C 
G(T)=RAMP L P(T)=SQRT(T) 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

INTEGER *4ITB( 12) /12«0/ 

REAL *4RTB(28) /28=«=0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(IOO), PR(IOO) 
ITB(3) = 5 
ITB(4) = 5 
WRITE (6,8) 

C 

DO 1 1=1,26 

READ (5,7) XX(I) , PR( I ) 

WRITE (6,7) XX(I),PR(I) 

1 CONTINUE 
C 

2 READ (5,7,EN0=4) TD,TFIN 
WRITE (6,5) TD,TFIN 
WRITE (6,6) 

X(l) = 0.0 
M = TFIN/TD 
A(l) = TD^*2 
A(2) = TO*=5=2 
A(3) = 0.0 
A(4) = 0.0 
C 

DO 3 N=1,M 
T = N*T0 
P = SQRT(T) 

FI = 2.0+P=4=TD 
F2 = -( 6.0 + P=4=TD) 

F3 = 6.0-P*TD 
F4 = (P*TD)-2.0 
C 

Q(N) = A(N)/F1 
C 

A(N + 1) = A(N+1)-(Q(N)=«=F2) 

A(N+2) = A(N+2)-( Q(N)«F3) 

A(N+3) = A(N+3)-(Q(N)*F4) 

A(N+4) = 0.0 
WRITE (6,7) T,Q(N) 

X(N) = T 

3 CONTINUE 
C 

ITB(l) = 1 
ITB(2) = 0 
RTB(6) = TD 
ITB(12) = 1 

CALL ORAWP ( 26,XX, PR, ITB, RTB) 

ITB(l) = 3 
IT8(2) = 2 

CALL ORAWP (M,X,Q,ITB,RTB) 

GO TO 2 

4 STOP 
C 

5 FORMAT Cl',' FOR DELTA T = • , F6 . 2 , 4X , ' TF I N= ' , F 4. 1 ) 

6 FORMAT (• t=',7X,' Q(N)=') 

7 FORMAT (2F10.5) 

8 FORMAT (• XX=',7X,» PR=') 

END 



PROGRAM 



8o 



dZl OIC ZIZ 020 




x-scale =1.0 units/inch 
y-scale = 0.5 units/inch 

FOR DELTA T = 0.2 

FIGURE ^19 
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1 



OTG 'JoG GCT 




FOR DELTA T = 0 



T 



0.3 

0.6 

0.9 

1.2 

1.5 
1.8 
2.1 
2.^ 
2.7 

3.0 
3.3 

3.6 
3.9 
H.2 

^.5 

i^.8 

5.1 



3 TPIN =5.2 

Q(N) 

0.0^158 

0.1551^ 

0.31700 

0 . 50^07 

0.69801 

0. 88639 

1.06213' 

1.22223 

1 . 366^11 

1. ^9595 
1.61280 
1.71905 

1.8l66i^ 

1.90718 

1.99200 

2.07212 

2.1iJ83^ 



FIGURE 50 
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'J2L, 010 




DELTA T = 0.4 


TFIN = 


T 


Q(N) 


0.4 


0.07102 


0.8 


0.25620 


1.2 


0.50068 


1.6 


0.75898 


2.0 


1.00322 


2.4 


1.22111 


2.8 


1.41066 


3.2 


1.57505 


3.6 


1.71916 


4.0 


1.84769 


4.4 


1.96446 


4.8 


2.07227 


5.2 


2.17312 



FIGURE 51 
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Q20 



-I 

O 



o 

— I 

o 



L> 

O 



c 

iT> 

O 




DELTA T = 0.5 


TFIN = 5.2 


T 


Q(N> 


0.5 


0.10622 


1.0 


0.36996 


1.5 


0.69092 


2.0 


1.00060 


2.5 


1.26987 


3.0 


1.49539 


3.5 


1.68484 


4.0 


1.84788 


4.5 


1.99233 


5.0 


2.12360 



FIGURE 52 
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OiO DIG 




FOR DELTA T = 0.6 TFIN =5.2 



T 



Q(N) 



0.6 

1.2 

1.8 

2.4 

3.0 

3.6 

4.2 

4.8 



0.14606 

0.49082 

0.87667 

1.21782 

1.49500 

1.71944 

1.90777 

2.07258 



FIGURE 53 
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ooti QiO Oi'j 020 




FOR DELTA T = 0.7 TFIN =5.2 



T 


Q(N). 


0.7 


0.18951 


1.4 


0.61452 


2.1 


1.05164 


2.8 


1.40840 


3.5 


1.68514 


4.2 


1.90808 


4.9 


2.09855 



FIGURE 54 
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I 






CZO 'jT.j OTO •'rJ', orc> 




FOR DELTA T = 0.8 TFIN = 5.2 



T 


Q(N) 


0.8 


0.23568 


1.6 


0.73797 


2.4 


1.21295 


3.2 


1.57461 


4.0 


1.84884 


4.8 


2.07304 



FIGURE 55 
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020 OTO OTO HOC 




FOR DELTA T = 0.9 



TFIN =5.2 



0.9 

1.8 

2.7 

3.6 

^.5 



Q(N) 

0.28383 

0.85903 

1.35987 

1.72042 

1.99350 



FIGURE 56 
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o'ic U2t; OiO Qi.'3 




FOR DELTA T = 1.0 TFIN =5,2 



T 


Q(N) 


1.0 


0.33333 


2.0 


0.97631 


3.0 


l.^<9299 


^.0 


1.85001 


5.0 


2. 12^161 



FIGURE 57 
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010 Oii 020 




FOR DELTA T = 1.2 TFIN =5.2 



T 


Q(N) 


1.2 


0.43^45 


2.H 


1.19662 


3.6 


1.72293 


H.8 


2.07450 



FIGURE 58 
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crjc: oiQ ois 




FOR DELTA T = 1.^1 TFIN =5.2 



T 


Q(N) 


1.4 


0.53603 


2.8 


1.39641 


4.2 


1.91440 



FIGURE 59 
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UOb OiO OiS 




FOR DELTA T = 1.8 TPIN =5.2 



T 



1.8 

3.6 



Q(N) 

0.73387 

1.73869 



FIGURE 60 
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0J.0 Oi'j 020 




FOR DELTA T = 2. TFIN = 5.2 



T 

2.0 



Q(N) 

0.82843 



FIGURE 6l 
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1 . Case 5 -Y+tY = t 

As case 2 follows case 1 , as calculations are 
concerned, case 5 follows case 4 with the only difference 
that P(t) Is now equal to t and this must be the value 
used to substitute in P before and during the division steps. 



Graphs, plots and solution values for several time 
Iterations follow (Pigs. 62 - 75 ). Computer algorithm for 
solution of this problem is shown In Program 5 (page 95 ) . 

FOR DELTA T = 0.2 TPIN = 5.2 



T 


Q(N) 




0.0 


0.0 


0.0 


0.2 


0.01961 


0.01980 


0.4 


0.07617 


0.07688 


0.6 


0.16332 


0.16473 


0.8 


0.27178 


0.27385 


1.0 


0.39094 


0.39347 


1.2 


0.51058 


0.51324 


1.4 


0.62220 


0. 62468 


1.6 


0.71991 


0.72186 


1.8 


0.80061 


0. 80209 


2.0 


0.86375 


0.86466 


2.2 


0. 91066 


0.91107 


2.4 


0.94380 


0.94386 


2.6 


0.96617 


0.96595 


2.8 


0.98040 


0.98016 


3.0 


0.98914 


0.98889 


3.2 


0.99424 


0.99402 


3.4 


0.99707 


0.99691 


3.6 


0.99857 


0.99847 


3.8 


0.99933 


0.99927 


4.0 


0.99970 


0.99966 


4.2 


0.99986 


0.99985 


4.4 


0.99993 


0.99994 


4.6 


0.99996 


0.99997 


4.8 


0.99997 


0.99999 


5.0 


0.99997 


1.00000 


Equation 


( 46 ) Is repeated 


here for convenience 




2-1 2-2 
T z + T z 




( 2 +PT) 


- ( 6 +PT)z~^ + ( 6 - 


PT)z“^ + (PT- 2 )z“^ 
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OUO' 



CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC3C 



Q'{T)=RAMP £ P(T)=T 



cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

INTEGER *4ITB( 12) /12*0/ 



c 


REAL *4RTB( 28)/28«0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(IOO), PR(IOO) 
ITB(3) = 5 
ITB(4) = 5 
WRITE (6,8) 


1 

c 

2 


DO 1 1=1,26 
READ (5,7) XX( I ) , PR( I ) 
WRITE (6,7) XX(I),PR(I) 
CONTINUE 

READ (5,7,END=4) TD,TFIN 
WRITE (6,5) TD,TFIN 
WRITE (6,6) 

X(l) = 0.0 
M = TFIN/TD 
A(l) = TD**2 
A(2) = TD*=!=2 
A(3) = 0.0 
A(4) = 0.0 


C 


DO 3 N=1,M 
T = N4=TD 
P = T 

FI = 2.0+P*TD 
F2 = -( 6.0+P*TD) 
F3 = 6.0-P*T0 
F4 = (P*TD)-2.0 


r 


Q(N) = A(N)/F1 


C 


A(N+1) = A(N+l)-( Q(N)*F2) 
A(N+2) = A(Ni-2)-(Q(N)=!=F3) 
A(N+3) = A(N^'3) -( Q(N)*F4) 
A(N+4) = 0.0 
WRITE (6,7) T,Q(N) 

X(N) = T 


3 

C 


CONTINUE 

ITB(l) = 1 
ITB(2) = 0 
RTB(6) = TD 
ITB(12) = 1 

CALL DRAWP ( 26 , XX , PR, I TB, RTB ) 
ITB( 1) = 3 
ITB(2) = 2 

CALL DRAWP ( M,X, Q, IT3, RTB) 

GO TO 2 



4 STOP 

5 FORMAT Cl',' FOR DELTA T = ' , F6. 2 , 4X , ' TF I N=' ,F4. 1 ) 

6 FORMAT (' T=',7X,' Q(N)=') 

7 FORMAT (2F10.5) 

8 FORMAT (' XX=',7X,' PR=') 

END 

PROGRAM 5 



95 



i 



PCD SCO iCO 




x-scale = 1.0 units/inch 
y-scale = 0.2 units/inch 

DELTA T = 0.2 



FIGURE 62 
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f 



PCO SCO 



-I — i — I — I — I — I — + 




FOR DELTA T = 0.3 TFIN = 5.2 



T 


Q(N) 


0.3 


0.0^1306 


0.6 


0.16158 


0.9 


0.32779 


1.2 


0.50723 


1.5 


0.67015 


1.8 


0.79871 


2.1 


0.88826 


2 . 1 \ 


0.9^372 


2.7 


0.97^36 


3.0 


0.989^8 


3.3 


0.99613 


3.6 


0.99873 


3.9 


0.99963 


^.2 


0.99990 


^.5 


0.99997 




0.99999 


5.1 


0.99999 



FIGURE 63 
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0^2 OO'i OQQ QO»i 




FOR DELTA T = 



T 



0. i] 

0.8 
1.2 
1.6 
2.0 

2.14 
2.8 

3.2 
3.6 
^4.0 

14.14 
H.8 

5.2 

FIGURE 



0.i| TFIN = 5.2 

Q(N) 

0.07^07 
0.2656^4 
0.50253 
0.71358 
0.86088 
0.9^360 
0.98120 
0.99^96 
0.99895 
0.999844 

0.99998 
1.00000 
1.00000 
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i 



i 



1 



fCA SOC '.'00 zoo 




DELTA T = 0.5 


TFIN = 


T 


Q(N) 


0.5 


0.11111 


1.0 


0.37776 


1.5 


0.66061 


2.0 


0.85859 


2.5 


0.956il9 


3.0 


0.99067 


3.5 


0.99876 


il.O 


0.99992 


^.5 


1.00000 


5.0 


1.00000 



FIGURE 65 
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I 



002 00 1 006 




0J2 


003 


FOR DELTA T = 


0 . 6 TFIN 


T 


Q(N) 


0.6 


0.1525^ 


1.2 


0.il8903 


1.8 


0.78765 


2.h 


0.9^321 


3.0 


0.99163 


3.6 


0.99960 


i|.2 


1.00001 


i|.8 


1.00000 



FIGURE 66 



QJ--. 

= 5.2 



100 



fOO SCO too ZT. 




FOR DELTA T = 0.7 TFIN =5.2 



T 


Q(N) 


0.7 


0.19679 


1.4 


0.59300 


2.1 


0.88036 


2.8 


0.98399 


3.5 


0.99986 


4.2 


1.00001 


4.9 


1.00000 



FIGURE 67 
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POO 500 ICO 200 




DELTA T = 0.8 


TFIN = 


T 


Q(N) 


0.8 


0.24242 


1.6 


0.68588 


2.4 


0.94231 


3.2 


0.99899 


4.0 


1.00011 


4.8 


0.99998 



FIGURE 68 
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iC'2 oo; vjs 




DELTA T = 0.9 


TFIN = 


T 


Q(N) 


0.9 


0.28826 


1.8 


0.76603 


2.7 


0.97993 


3.6 


1.00165 


4.5 


0.99966 



FIGURE 69 
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PCO 600 IGC ?0C 




FOR DELTA T = 1.0 TFIN = 5.2 



T 


Q(N) 


1.0 


0.33333 


2.0 


0.83333 


3.0 


1.00000 


il.O 


1.00000 


5.0 


1.00000 



FIGURE 70 



10^1 






I 



POO 900 too rcTi 




FOR DELTA T = 1.2 TFIN =5.2 



T 



Q(N) 



1.2 

2.H 

3.6 

n.Q 



0. ill 860 
0.93328 

1.00929 

0.99722 



FIGURE 71 
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i 



9C0 !'C0 aoo 




FOR DELTA T = 1 . ij TFIN = 5.2 



T 



Q(N) 



1.4 

2.8 

4.2 



0.49495 

0.99659 

1.00083 



FIGURE 72 
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HH HH HH HH HH 

HH HH HH HH HH 

; »» »» »» »» »» 

. 



CuO 900 KO zy\ 







FOR DELTA T = 1.6 TFIN =5.2 



T 


Q(N) 


1.6 


0.56140 


3.2 


1.03450 


4.8 


0.98888 



FIGURE 73 
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POO SCO !CO SCO 







FOR DELTA T = 1.8 TFIN =5.2 



T 

1.8 

3.6 



Q(N) 

0.61832 

1.05581 



FIGURE 7^ 
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foo SCO 



+ 




FOR DELTA T = 2.0 TFIN =5.2 



T Q(N) 



2.0 0.66667 

4.0 1.06667 



FIGURE 75 
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2 . 



Case 6 - Y + t^Y = t 

Calculations follow exactly Case 4 except for ■ 

2 

P(t) = t . In the division steps, substitution of P for 
2 

its value t must be done in this case. 

Solutions for several time iterations follow 
(Figs. 76 to 89). Computer solution is in Program 6 (page 
111 ). 



FOR DELTA T = 0.2 TPIN - 5-2 



T 


Q(N) 


PR 


0.2 


0.01992 


0.01997 


0.4 


0.07858 


0.07898 


0.6 


0.17116 


0.17243 


0.8 


0 . 2 8666 


0.28923 


1.0 


0.40755 


0.41144 


1.2 


0.51294 


0.51742 


1.4 


0.58451 


0.58837 


1.6 


0.61301 


O.6I523 


1.8 


0.60127 


0.60163 


2.0 


0.56176 


0.56095 


2.2 


0.51015 


0.50914 


2.4 


0.45891 


0.45829 


2.6 


0.41443 


0.41425 


2.8 


0.37796 


0.37801 


3.0 


0.34824 


0. 34833 


3.2 


0.32354 


0. 32360 


3.4 


0.30253 


0.30257 


3.6 


0.28433 


0.28436 


3.8 


0.26835 


0.26837 


4.0 


0.25418 


0.25420 


4.2 


0.24150 


0.24152 


4.4 


0.23008 


0.23009 


4.6 


0.21973 


0.21974 


4.8 


0.21029 


0.21030 


5.0 


0.20166 


0.20167 


5.2 


0.19372 





Equation ( 46 ) is repeated here for convenience: 
*(z) = 



T^z~^ + T^z~^ 



(2 + PT) - ( 6 +PT)z ^+( 6 -PT)z ^ + (PT- 2 )z ^ 
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CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC3C 

c 

C C(T)=RAMP £ P(T)=T**2 

C 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

INTEGER =i'4ITB( 12)/12*0/ 

REAL *4RTB(28)/28*0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(IOO), PR(IOO) 
ITB(3) = 5 
ITB(4) = 5 
WRITE (6,8) 

C 

DO 1 1=1,26 

READ (5,7) XX( I) ,PR(I) 

WRITE (6,7) XX(I),PR(I) 

I CONTINUE 



READ (5,7,END=^) TD,TFIN 
WRITE (6,5) TD,TFIN 
WRITE (6,6) 

X(l) = 0.0 
M = TFIN/TD 



A(l) 

A(2) 

A(3) 

A(4) 



TD**2 

TD--2 

0.0 

0.0 



DO 3 N=1,M 
T = N^f'TD 

p = 

FI = 2.0+P*TD 
F2 = -( 6.0+P*TD) 

F3 = 6. 0-P*TD 
F4 = (P=!=TD)-2.0 

Q(N) = A(N)/F1 

A(N+1) = A(Ntl)-( Q(N)*F2) 
A(N+2) = A( N+2) -(Q(N)*F3) 
A(N+3) = A(iN+3) -(Q(N)«F4) 
A(N+4) = 0.0 
WRITE (6,7) T,Q(N) 

X(N) = T 
CONTINUE 



ITB(l) = 1 
ITB(2) = 0 
RTB(6) = TD 
ITB(12) = 1 
CALL DRAWP 
IT3(1 ) =3 
ITB(2) = 2 
CALL DRAWP 
GO TO 2 
STOP 



( 26,XX,PR,ITB,RTB) 



( M,X,Q ,ITB,RTB) 



5 

6 
7 



FORMAT 
FORMAT 
FORMAT 
8 FORMAT 
END 



('!»,' FOR 
(' T=',7X, 

( 2F10.5) 

(' XX=',7X 



DELTA T=’ ,F6.2,4X,'TFIN=* ,F4.1) 
' Q(N)=’) 

,' PR=M 



PROGRAM 6 
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x-scale = 1.0 iinlts/inch 
y-scale =0.2 units/inch 

FOR DELTA T = 0.2 



FIGURE 76 
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t 



i-cz oa; oc'6 crjJ 







FOR DELTA T = 0.3 



T 



0.3 

0.6 

0.9 

1.2 

1.5 
1.8 
2.1 

2.4 

2.7 

3.0 
■3.3 

3.6 
3.9 
4.2 

4.5 

4.8 

5.1 



TFIN =5.2 

Q(N) 

0.04440 

0.16964 

0.34372 

0.50736 

0.60020 

0.60070 

0.53792 

0.45981 

0.39529 

0.34810 

0.31257 

0.28430 

0.26104 

0.24149 

0.22477 

0.21028 

0.19760 



FIGURE 77 
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0C2 OQ'i OOS OOU 




FOR DELTA T = 0 



T 



0.4 

0.8 

1.2 

1.6 

2.0 

2.4 

2.8 

3.2 
3.6 
4.0 

4.4 
4.8 

5.2 



4 TPIN =5.2 

Q(N) 

0.07752 

0.27929 

0.49964 

0.60565 

0.56420 

0.46136 

0.37768 

0.32332 

0.28429 

0.25415 

0.23007 

0.21028 

0.19371 



FIGURE 78 



Il4 



coo SQO too ZCT. 




DELTA T = 0.5 


TFIN = 


T 


Q(N) 


0.5 


0.11765 


1.0 


0.38824 


1.5 


0.58635 


2.0 


0.56576 


3.0 


0.34709 


3.5 


0.29320 


4.0 


0.25405 


4.5 


0.22480 


5.0 


0.20162 



FIGURE 79 
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0C2 001 006 003 




FOR DELTA T = 0.6 TFIN =5.2 



T 


Q(N). 


0.6 


0.16245 


1.2 


0.47829 


1. 8 


0.59415 


2.4 


0.46797 


3.0 


0.34576 


3.6 


0.28482 


4.2 


0.24117 


4.8 


0.21043 



FIGURE 80 
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I 



ICZ mi 0‘J5 OQt? 




FOR DELTA T = 0.7 TFIN = 5.2 



T 


Q(N) 


0.7 


0.20913 


l.il 


0.53871 


2.1 


0.5^812 


2.8 


0.37850 


3.5 


0.29218 


^.2 


0 . 2 ijl 77 


^.9 


0.20568 



FIGURE 8l 
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POO 000 too 200 




FOR DELTA T 



0.8 TFIN =5.2 



T 



Q(N) 



0.8 
1.6 
2 A 
3.2 

l\.S 



0.25^78 
0.56796 
0. il801iJ 

0.31670 

0.25669 

0.20888 



FIGURE 82 
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0C2 TOi 006 003 




FOR DELTA T = 0.9 TPIN = 5.2 



T 


Q(N) 


0.9 


0.29681 


1. 8 


0.57104 


2.7 


0.41198 


3.6 


0.27744 


4.5 


0.22788 



FIGURE 83 
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i 

's 

I 



E 




, DELTA T = 1.0 


TFIN = 


T 


Q(N) 


1.0 


0.33333 


2.0 


0.55556 


3.0 


0.3535^ 


i|.0 


0.251^0 


5.0 


0.20298 



FIGURE 
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pro sot\ '-00 ZcA 







FOR DELTA T = 1.2 TFIN = 5.2 
T Q(N) 

1.2 0.38627 

2.4 0.49653 

3.6 0.27125 

4.8 0.21600 



FIGURE 85 
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I 



0C2 001 006 




I DELTA T = 1.4 


TPIN = 5.2 


T 


Q(N) 


1.4 


0.41315 


2.8 


0.42946 


4.2 


0.22270 


FIGURE 86 
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LC2 000 OOG OOi? 




FOR DELTA T = 



T 

1.6 

3.2 

^.8 



1.6 TFIN = 5.2 

Q(N) 

0. ill995 

0.36988 

0.192il6 



FIGURE 87 
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i 



PCO &CO ICO 7.Z0 coc 




FOR DELTA T = 1.8 TFIN = 5.2 

T Q(N) 

1.8 0.ill369 

3.6 0.32118 



FIGURE 88 



124 



y^2 QQ'l 0^6 OQJ 




FOR DELTA T = 2.0 TFIN = 5.2 

T Q(N) 

2.0 0.40000 

4.0 0.28235 



FIGURE 89 
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' Dii' I"; ■"! 

i ::: 

“■! IS! IS! IS! IS! IS! : IS; 




D. CASE 7-Y + \/FY = t^ 



Following the procedure as in II C : 

II + PY = P 

The Laplace transform in powers of s“^: 

sy + Py = 2/s^ y = — — 

1 + Ps"-^ 



Substituting the z-forms of Table I into the expression 
above yields : 



y^*(z) 



z”^ + + z~^ 

(1-z-l)^ 



T + P I (1+ Z)~^ 
^ ^ 2 • ^-l 



(1 - z) 



T^ 

720 



Dividing by T and rearranging in crescent pov^ers of z 



- 1 , 



(z) = 



-T^ + (12ilT3)z ^ + (47^T^)z ^ + 12^T^z ^ = T^z~^ 



360+l80PT-(li]iJ0+360Pr)z ^+2l60z ^-(li4i40-360PT)z 



-34 



(i*7) 



+ (360-180PT)z 






Choosing for first value for T = 0.2, and substituting 
in the previous equation, the denominator is divided into 
the numerator, after the value of P is replaced in the 
expression above for the corresponding time. 
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The results and the plots for the precise and approximated 
solution follow, with iteration times from 0.2 sec. to 
2.0 sec. (Pigs. 90 to 103 ). Computer algorithm for solution 
of this problem is in Program 7 (page 128 ). 



DELTA T = 0.2 TFIN = 5.2 



T 


Q(N) 




0.0 


- 0.00002 


0.0 


0.2 


0.00255 


0.00261 


0.4 


0.01985 


0.02018 


0.6 


0.06428 


0.06511 


0.8 


0.14500 


0.14648 


1.0 


0.26791 


0.27012 


1.2 


0.43602 


0.43897 


1.4 


0.65001 


0.65363 


1.6 


0.90881 


0.91298 


1.8 


1.21016 


1.21470 


2.0 


1.55110 


1.55590 


2.2 


1.92838 


1.93330 


2.4 


2.33869 


2.34370 


2.6 


2.77893 


2.78380 


2.8 


3.24627 


3.25100 


3.0 


3.73821 


3.74280 


3.2 


4.25262 


4.25700 


3 .^ 


4.78767 


4.79180 


3.6 


5.34186 


5.34590 


3.8 


5.91392 


5.91780 


4.0 


6.50283 


6.50670 


4.2 


7.10775 


7.11160 


4.4 


7.72795 


7.73200 


4.6 


8.36285 


8.36720 


4.8 


9.01195 


9.01680 


5.0 


9.67481 


9.68030 



127 



cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

C Q(T) = J**Z & P(T) =SQRT(T) 

c 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

INTEGER *4ITB( 12)/12=^=0/ 

REAL *4RTB(28)/28*0.0/ 

DIMENSION A(IOO), QllOO), XIIOO), XX(IOO), PR(IOO) 
ITBI3) = 5 
ITBI4) = 5 
WRITE (6,8) 

C 

DO 1 1=1,26 

READ (5,5) XX(I),PR(I) 

WRITE (6,5) XX(I),PR(I) 

1 CONTINUE 
C 

2 READ (5,5,END=4) TD,TFIN 
WRITE (6,6) TD,TFIN 
WRITE (6,7) 

M = TFIN/TD 
T3 = TD*^-3 
A(l) = -T3 
A(2) = 124.0*T3 
A(3) = 474.0=!'T3 
A(4) = 124.0*T3 
A(5) = -T3 
C 

DO 3 N=1,M 
T = (N-l)'^TD 
P = SQRT( T) 

PT = P*TD 

FI = 360. 0+180. 0*PT 
F2 = -( 1440. 0+360. 0*PT) 

F3 = 2160.0 

F4 = -{ 1440. 0-360. 0-PT) 

F5 = 360. 0-180. 0=t=PT 



Q(N) = A(N)/F1 

A(N+1) = A(N + 1)-(Q(N)=!'F2) 
A(N+2) = A( N+2) -( Q(N)*F3) 
A(N+3) = A(N+3)-(Q(N)*F4) 
A(N+4) = A(N+4)-(Q(N)*F5) 
A(N+5) = 0.0 
WRITE (6,5) T,Q(N) 

X(N) = T 
CONTINUE 



ITB(i) = 1 
ITB(2) = 0 
ITB(12) = 1 
CALL DRAWP 
ITB( 1 ) = 3 
ITB(2) = 2 
CALL DRAWP 
GO TO 2 
STOP 



(26,XX,PR,ITB,RTB) 



( H,X, Q, I TB, RTB) 



FORMAT 

FORMAT 

FORMAT 

FORMAT 

END 



( 2F10.5) 

{'I', 'FOR DELTA T= ’ , F6 . 2 , 4X , ’ TF I N= ’ , F4. 1 ) 
{ • T= • ,8X, • 0(N)=* ) 

(• XX=',7X,‘ PR=') 



PROGRAM 7 
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fiCO SCO I CO 200 




x-scale =1.0 units/inch 
y-scale = 2.0 units/inch 

FOR DELTA T = 0.2 



FIGURE 90 
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i 



I 



I 

i 



UC2 Q'.'-i 006 



g 




FOR DELTA T = 0 



T 

0.0 

0.3 

0.6 

0.9 

1.2 

1.5 
1.8 
2.1 
2.i| 
2.7 

3.0 

3.3 

3.6 
3.9 
H.2 
^.5 



3 TFIN =5.0 
Q(N) 

-0.00007 

0.00832 

0.06328 

0.19866 

0.43235 

0.76902 

1.20437 

1.72913 

2.33229 

3.00312 

3.73225 

4.51208 

5.33673 

6.20176 

7.10386 

8.04059 



FIGURE 91 
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f-CO SCO vco z:r. 




FOR DELTA T = 0.4 


TFIN = 


T 


Q(N) 


0.0 


- 0.00018 


0.4 


0.01894 


0.8 


0.14071 


1.2 


0.42726 


1.6 


0.89625 


2.0 


1.53631 


2.4 


2.32327 


2.8 


3.23134 


3.2 


4.23875 


3.6 


5.32927 


4.0 


6.49164 


4.4 


7.71831 



FIGURE 92 
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i 



* 



FCO 6C0 




DELTA T = 0.5 


TFIN =5.0 


T 


Q(N) 


0.0 


-0.00035 


0.5 


0.03541 


1.0 


0.25665 


1.5 


0.75346 


2.0 


1.52519 


2.5 


2.52824 


3.0 


3.71282 


3.5 


5.03967 


4.0 


6.48280 


4.5 


8.02638 



FIGURE 93 
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' f 



n 

! ' 
h 

II 

i 



OQC UC2 00 ; 005 00*3 




FOR DELTA T = 0.6 TFIN =5.0 
T Q(N) 



0.0 


-0.00060 


0.6 


0.05842 


1.2 


0.41309 


1.8 


1.17315 


2. il 


2.29728 


3.0 


3.69936 


3.6 


5.30783 


4.2 


7.07854 



FIGURE 9'^ 
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PCO &CO 




FOR DELTA T = 0.7 TPIN = 5.0 



T 



Q(N) 



0.0 

0.7 

1.4 
2.1 
2.8 

3.5 
4.2 



-0.00095 

0.08844 

0.61018 

1.67849 

3.18952 

5.01324 

7.06634 



FIGURE 95 
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i! 



I 



0C2 00'. QOS 003 




FOR DELTA T = 0.8 TFIN = 5.0 



T 



Q(N) 



0.0 

0.8 

1.6 

2 .^ 

3.2 

4.0 



-0.00142 
0.12570 
0.84689 
2 .26036 
4.18196 
6.44458 



FIGURE 96 
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PCO &CO 




FOR DELTA T 
T 

0.0 

0.9 

1.8 

2.7 

3.6 



=0.9 TFIN =5.0 
Q(N) 

- 0.00202 

0.17030 

1.12150 

2.90981 

5.25900 



FIGURE 97 
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!JC2 QCl'l OGS ooa 




FOR DELTA T 
T 

0.0 

1.0 

2.0 

3.0 

.0 



= 1.0 TFIN = 5.0 
Q(N) 

- 0.00278 

0.22222 

l.i»3192 

3.61877 

6.iJ0907 



FIGURE 98 
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oo<3 ooa 




FOR DELTA T = 1.2 TFIN =5.0 



T Q(N) 



0.0 


-0.00480 


1.2 


0.34756 


2.4 


2.15138 


3.6 


5.18870 



FIGURE 99 
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I 



POO 900 VCO ZZO 




FOR DELTA T = 1.4 TFIN =5.0 



T 



Q(N) 



0.0 

1.4 

2.8 



-0.00762 

0.50030 

2.98831 



FIGURE 100 
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PCO 9Q0 ^CO 




FOR DELTA T = 1.6 


TFIN = 


T 


Q(N) 


0.0 


-0.01138 


1.6 


0.67862 


3.2 


3.92797 



FIGURE 101 



1^»0 






I 



'JC2 QOl QOS 00^ 




FOR DELTA T = 1.8 TPIN = 5.0 

T Q(N) 

0.0 - 0.01620 

1.8 0.88064 



FIGURE 102 



I4l 



1 










+ 



o:; 0J2 QC3 001 

FOR DELTA T = 2.0 TFIN =5.0 



T 


Q(N) 


0.0 

2.0 


-0.02222 

1.10457 



FIGURE 103 
li|2 



1. Case 8 - Y + tY = 

Same calculations as in case 7 with exception for 
the substitution on the value of P during the division 
steps, of Equation (^7), repeated here, which now has the 
value P = t . 

The algorithm for computer solution (Program 8, 
page l^l^l), the output for several time Iterations and the 
respective graphs follow (Pigs. 10^1 to 117) 



-T^+12^T^)z V(^7^T^z ^+12 ^-T^z ^ 

360+l80PT-(lilil0+360PT)z“^+2l60z"^-(lilil0-360PT)z“^. . . 

•• •+(360 - 180 FT)z"^ 



FOR DELTA T = 0.2 TFIN = 5.0 



T 


Q(N) 


m 


0.0 


-0.00002 


0.0 


0.2 


0.00261 


0.00264 


0. 


0.02041 


0.02066 


0.6 


0.06629 


0.06707 


0.8 


0.14905 


0.15068 


1.0 


0.27254 


0.27522 


1.2 


0.43567 


0.43942 


l.l\ 


0.63339 


0.63801 


1.6 


0.85809 


0.86328 


1.8 


1.10124 


1.10660 


2.0 


1.35473 


1.35990 


2.2 


1.61184 


1.61660 


2. A 


1.86765 


1. 87190 


2.6 


2.11908 


2.12280 


2.8 


2.36463 


2.36780 


3.0 


2.60393 


2.60670 


3.2 


2.83735 


2.83980 


3.^ 


3.06563 


3.06780 


3.6 


3.28958 


3.29160 


3.8 


3. 51003 


3.51190 


. 0 


3.72765 


3. 72940 


^.2 


3.94301 


3.94460 


4.4 


4 .15655 


4.15810 


4.6 


4.36861 


4.37000 


4.8 


4.57945 


4.58080 



1^13 



non 



cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

Q(T) = T**2 £1 P(T) =T 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

INTEGER «4ITB( 12)/12«0/ 

REAL *4RTB(28)/28«0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(IOO), PR(IOO) 
ITB(3) = 5 
ITB(4) = 5 
WRITE (6,8) 

C 

DO 1 1=1,26 

READ (5,5) XX(I),PR(I) 

WRITE (6,5) XX(I),PR(I) 

1 CONTINUE 
C 

2 READ (5,5,END=4) TD,TFIN 
WRITE (6,6) TD,TFIN 
WRITE (6,7) 

M = TFIN/TD 
T3 = JD^^3 
A(l) = -T3 
A(2) = 124.0=s=T3 
A(3) = 474.0=5=T3 
A(4) = 124.0*T3 
A(5) = -T3 
C 

DO 5 N=1,M 
T = (N-l)-TD 
P = T 
PT = P^TD 

FI = 360. 0 + 180. 0*PT 
F2 = -( 1440. 0+360. 0*PT) 

F3 = 2160.0 

F4 = -( 1440. 0-360. 0*PT) 

F5 = 360. 0-180. 0*PT 
C 

Q(N) = A(N)/F1 
C 

A(N+1) = A(N+1) -(Q(N)*F2) 

A(N+2) = A( N+2)-( Q(N)*F3) 

A(N+3) = A(N+3)-(Q(N)*F4) 

A(N+4) = A(N+4)-(Q(N)*F5) 

A(N+5) = 0.0 
WRITE (6,5) T,Q(N) 

X(N) = T 

3 CONTINUE 
C 

ITB(l) = 1 
ITB(2) = 0 
ITB(12) = 1 

CALL DRAWP ( 26 , XX , PR , I Tb,RTB ) 

ITB( 1) = 3 
ITB(2) = 2 

CALL DRAWP ( M ,X , Q , I TB , RTB ) 

GO TO 2 

4 STOP 
C 

5 FORMAT (2F10.5) 

6 FORMAT Cl', 'FOR DELTA T= ’ , F6 . 2 , 4X , ' TF I N = * , F4 . 1 ) 

7 FORMAT (' T=',8X,' 0(N)=') 

8 FORMAT C XX=',7X,' PR=') 

END 



PROGRAM 8 



20D TOO 




x-scale = 1.0 units/inch 
y-scale =1.0 units/inch 

FOR DELTA T = 0.2 

FIGURE 10i< 






0Q2 




FOR DELTA T = 0.3 



TFIN = 5.0 



Q(N) 



0.0 


- 0.00007 


0.3 


0.00861 


0.6 


0.06534 


0.9 


0.20305 


1.2 


0.43105 


1.5 


0.73670 


1. 8 


1.09444 


2.1 


1.47676 


2.4 


1.86218 


2.7 


2.23828 


3.0 


2.60051 


3.3 


2.94932 


3.6 


3.28723 


3.9 


3.61706 


4.2 


3.94109 


4.5 


4.26094 



FIGURE 105 



1^16 



2C0 




FOR DELTA T = 0.^1 
T 

0.0 

0 . 

0.8 

1.2 

1.6 

2.0 

2 .^ 

2.8 

3.2 

3.6 

i ].0 



TFIN =5.0 
Q(N) 

- 0.00018 
0.01975 
0. li]i|40 
0. il 2 il 70 
0.84250 
1.33874 
1.85444 

2.35490 

2.83016 

3.28380 

3.72261 

4.15198 



FIGURE 106 
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I 



- 






ZuO Tcn cno 




FOR DELTA T = 0.5 


TFIN = 


T 


Q(N) 


0.0 


- 0.00035 


0.5 


0 . 0370 ^ 


1.0 


0.25926 


1.5 


0.71717 


2.0 


1.32660 


2.5 


1.97228 


3.0 


2.58929 


3.5 


3.167^^ 


i|.0 


3.71877 


^.5 


i]. 25^87 



FIGURE 107 












I 



too 200 TOO COO 




FOR DELTA T = 0.6 TFIN =5.0 



T Q(N) 



0.0 


-0.00060 


0 . 6 


0.06102 


1.2 


0.i^0738 


1.8 


1.05761 


2.M 


1.83168 


3.0 


2.58151 


3.6 


3.27^10 


A. 2 


3.93011 



FIGURE 108 






zno Tcc 




FOR DELTA T = 0.7 TFIN =5.0 



T Q(N) 



0.0 


-0.00095 


0.7 


0.09183 




0.58367 


2.1 


1. i^2363 


2.8 


2.32707 


3.5 


3.155^4 


i ^.2 


3.92^82 



FIGURE 109 
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FOR DELTA T = 0.8 



TFIN =5.0 



T 



Q(N) 



0.0 

0.8 

1.6 

2.4 

3.2 

4.0 



-0.00142 

0.12929 

0.78206 

1.79807 

2.80113 

3.70243 



FIGURE 110 
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CCO 200 TOC 000 




022 


uuo 


FOR DELTA T = 0.9 


TFIN = 5 


T 


Q(N) 


0.0 


=0.00202 


0.9 


0.17295 


1.8 


0.99663 


2.7 


2.16991 


3.6 


3.25361 • 



FIGURE 111 
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00*: oci ooi? 0C3 



o 

o 




DELTA T = 1.0 


TFIN = 1 


T 


Q(N) 


0.0 


-0.00278 


1.0 


0.22222 


2.0 


1.22222 


3.0 


2.53333 


i^.O 


3.68888 



FIGURE 112 
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ex zee TCO 




FOR DELTA T = 1.2 TFIN = 5.0 
T Q(N) 



0.0 


-0.00480 


1.2 


0.33488 


2.4 


1.69089 


3.6 


3.22784 



FIGURE 113 
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QQ'l 




FOR DELTA T = 1 . TPIN = 5.0 
T Q(N) 



0.0 


- 0.00762 


1.4 


0.46195 


2.8 


2.16618 



FIGURE 114 
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OL'O QQi. 0Q2 UC3 




FOR DELTA T = 1.6 TPIN =5.0 



T 



Q(N) 



0.0 

1.6 

3.2 



- 0.01138 

0.59883 

2.63754 



FIGURE 115 
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TOO ?00 TCO XO 




FOR DELTA T = 1.8 TPIN =5.0 

T Q(N) 

0.0 - 0.01620 

1.8 0 . 7^198 



FIGURE Il6 



157 



OUC QOl 002 U03 001 




FOR DELTA T = 2.0 TFIN =5.0 

T Q(N) 

0.0 - 0.02222 

2.0 0.88889 



FIGURE 117 
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2 . Case 9 — Y + t^Y = 

Same reasoning as in case 8 . This time P = t^. 
Computer algorithm for solution of this problem is 
shown in Program 9 (Page i 6 o) . Solutions and graphs for 
several T’s follow in Figures 118 to I3I. 



FOR DELTA T = 0.2 TFIN = 5.0 



T 


Q(N) 


PR 


0.0 


-0.00002 


0.0 


0.2 


0.00266 


0.00266 


0.4 


0.02098 


0.02110 


0.6 


0.06883 


0.06947 


0.8 


0.15509 


0.15690 


1.0 


0.27985 


0.28347 


1.2 


0.43228 


0.43785 


1.4 


0.59256 


0.59934 


1.6 


0.73812 


0.74469 


1.8 


0.85183 


0.85686 


2.0 


0.92749 


0.93051 


2.2 


0.96977 


0.97125 


2.4 


0.98924 


0.99002 


2.6 


0.99645 


0.99714 


2.8 


0.99857 


0.99934 


3.0 


0.99909 


0.99988 


3.2 


0.99924 


0.99998 


3.4 


0.99932 


1.00000 


3.6 


0.99938 


1.00000 


3.8 


0.99943 


1.00000 


4.0 


0.99947 


1.00000 


4.2 


0.99950 


1.00000 


4.4 


0.99952 


1.00000 


4.6 


0.99954 


1.00000 


4.8 


0.99956 


1.00000 



Equation ( 47 ) is repeated for convenience: 



-T^+ ( 124T^) z~^+ ( 474T^) z~^+124t^z~^-T^z~^ 

360 +l 80 FT-(l 440 + 360 PT)z“^+ 2 l 60 z“^-(l 440 - 360 PT)z“ 3 +- • • 

...+( 360-1 80 fT)z"^ 
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cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

C Q(T) = T**2 G P(T) = 1**2 
C 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

INTEGER *4ITB( 12)/12*0/ 

REAL *4RTB( 28 )/28*0.0/ 

DIMENSION A(IOO), Q(IOO), X(IOO), XX(IOO), PR(IOO) 
ITB(3) = 5 
ITB(4) = 5 
WRITE (6,8) 

C 

DO 1 1=1,26 

READ (5,5) XX(I) ,PR(I) 

WRITE (6,5) XX(I),PR(I) 

1 CONTINUE 
C 

2 READ (5,5,END=4) TD,TFIN 
WRITE (6,6) TD,TFIN 
WRITE (6,7) 

M = TFIN/TD 
T3 = TD^^^'S 
■ A(l) = -T3 

A(2) = 124.0*T3 
A(3) = 474.0*T3 
A(4) = 124.0*T3 
A(5) = -T3 
C 

DO 3 N=1,M 
T = (N-1)*TD 
P = T=!=*2 
PT = P*TD 

FI = 360. 0 + 180. 0=«=PT 
F2 = -( 1440. 0+360. 0*PT) 

F3 = 2160.0 

F4 = -( 1440. 0-360. 0*PT) 

F5 = 360. 0-180. 0*PT 
C 

Q(N) = A(N)/F1 
C 

A(N+1) = A(N+1)-(Q(N)*F2) 

A(N+2) = A(N+2)-(Q(N)-F3) 

A(N+3) = A(N+3)-( Q(N)-F4) 

A(N+4) = A( N+4) -(Q(N)*F5) 

A(N+5) = 0.0 
WRITE (6,5) T,0(N) 

X(N) = T 

3 CONTINUE 
C 

ITB(l) = 1 
ITB(2) = 0 
ITB(12) = 1 

CALL DRAWP ( 26 , XX , PR , I TB, RT B ) 

ITB(l) = 3 
ITB(2) = 2 

CALL DRAWP ( M , X , Q , I TB , RTB ) 

GO TO 2 

4 STOP 
C 

5 FORMAT (2F10.5) 

6 FORMAT ('1','FOR DELTA T= • , F6 . 2 , 4X , ' TF I N= ' , F4. 1 ) 

7 FORMAT (' T=',8X,' Q(N)=') 

8 FORMAT (• XX=',7X,* PR=') 

END 



PROGRAM 9 



OQ'l 006 006 




x-scale = 1.0 units/inch 
y-scale =0.2 unlts/lnch 

FOR DELTA T = 0.2 



FIGURE 118 
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OOC u:2 03'". 036 009 



^ 1 1 1 h — I h 



Uuu 



oa 2 


C'u.j 


DELTA T = 0 . 3 


TFIN = 


T 


Q(N) 


0.0 


- 0.00007 


0.3 


0.00888 


0.6 


0.06808 


0.9 


0.20990 


1.2 


0.42549 


1.5 


0.65988 


1.8 


0.84534 


2.1 


0.94944 


2.4 


0.98840 


2.7 


0.99809 


3.0 


0.99819 


3 . 3 


0.99851 


3.6 


0.99877 


3.9 


0.99895 


4.2 


0.99910 


4.5 


0.99921 



oa-j 



FIGURE 119 
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eoo 900 VuO 20 C 30C 




FOR DELTA T = 0.4 TPIN =5.0 



T 



Q(N) 



0.0 

0.4 

0.8 

1.2 

1.6 

2.0 

2.4 
2.8 
3.2 
3.6 
4.0 

4.4 



- 0.00018 

0.02067 

0.15013 

0.41634 

0.71810 

0.91765 

0.98739 

0.99659 

0.99714 

0.99787 

0.99827 

0.99859 



FIGURE 120 
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PCO 900 t-CO 2CC 




: DELTA T = 0.5 


TFIN = ' 


■T 


Q(N) 


0.0 


-0.00035 


0.5 


0.03922 


1.0 


0.26275 


1.5 


0.63278 


2.0 


0.90925 


2.5 


0.99187 


3.0 


0.99^99 


3.5 


0.99641 


i^.O 


0.99731 


^.5 


0.99789 



FIGURE 121 
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Q'Jl OOG OQ'J 




: DELTA T = 0.6 


TFIN = ! 


T 


Q(N) 


0.0 


-0.00060 


0.6 


0.06498 


1.2 


0.39243 


1. 8 


0.80675 


2.4 


0.98482 


3.0 


0.99326 


3.6 


0.99498 


4.2 


0.99658 



FIGURE 122 
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0C2 OQl 006 



I- 




FOR DELTA T = 0.7 TFIN =5.0 



T Q(N) 



0.0 


-0.00095 


0.7 


0.09760 


1.4 


0.52265 


2.1 


0.91859 


2.8 


0.99655 


3.5 


0.99033 


4.2 


0.99646 



FIGURE 123 
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001 006 00 J 




DELTA T = 0.8 


TPIN =5.0 


T 


Q(N) 


0.0 


-0.00142 


0.8 


0.13588 


1.6 


0.64020 


2 A 


0.97679 


3.2 


0.98920 




0.99299 



FIGURE 12^1 
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i 



PCO 000 ICO 




FOR DELTA T = 0.9 TFIN = 5-0 



T 


Q(N) 


0.0 


-0.00202 


0.9 


0.17809 


1.8 


0.73807 


2.7 


0.99964 


3.6 


0.98234 



FIGURE 125 
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I 



'I 



UC 2 Q!j<a OL'9 ooa 




DELTA T = 1.0 


TFIN = : 


T 


Q(N) 


0.0 


-0.00278 


1.0 


0.22222 


2.0 


0.81481 


3.0 


1.00337 


4.0 


0.98017 



FIGURE 126 
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PCO GOO too 200 




: DELTA T = 1.2 


TFIN = 


T 


Q(N) 


0.0 


-0.00480 


1.2 


0.30901 


2.4 


0.91428 


3.6 


0,99117 



FIGURE 127 
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i 



00 - 



0:i Q02 


SOS 00- 


FOR DELTA T = 1.^1 


TFIN =5.0 


T 


Q(N) 


0.0 
1. H 
2.8 


-0.00762 

0.38561 

0 . 96 ^ 17 ^ 


FIGURE 128 
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f'OO SCO I'CO zco 




FOR DELTA T = 1.6 TPIN = 5.0 



T Q(N) 



0.0 


-0.01138 


1.6 


0. 


3.2 


0.98867 



FIGURE 129 
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POO OCO 100 200 




FOR DELTA T = 1.8 TFIN =5.0 



T 



0.0 

1.8 



Q(N) 

- 0.01620 

0.49642 



FIGURE 130 
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JC2 001 005 008 




FOR DELTA T = 2.0 TFIN =5.0 

T Q(N) 

0.0 -0.02222 

2.0 0.53333 



FIGURE 131 
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E. CASE 10 



No one of the previous cases has a time varying 
coefficient which Is exponential. So, from Ref. 2 this 
case was extracted. The equation Is: 

C || + = 0 being r(t) = R exp (i^) 

where Is a factor In the variation of r. Dividing all 
the terms by C yields: 



de 

dt 



+ 



1 f-t y. 

FfT- exp(f^) 
±2 ■ !]_ 



e = 0 



where T 2 = RC . 

For this example and T 2 were made equal to one and 
= 1.0 . Rewriting the equation and expressing its 
Laplace transform in powers of s”^: 



e + Ce = 0 C = e"'^ 

sE - 1 + CE = 0 

F = 1 = i _ Cs~^ 

® ^ ® 1 + Cs"^ 

If the rules In Example 2, Section IIB v;ere followed, 
substituting the z-forms of s“^, dividing by T and adding 
the final results of the quotient to 1.0 which Is the time 
function for 1/s: 
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®A* 



^ 1 + lOz”^ + z"^ 



1 + C 



T 1 + z”^ 



For T = 0.1 



F * 



-0.1- z“^ - 0. Iz”^ 



(12/C + 0.6) - (24/C)z"^ + (12/C - 0.6)z 



-2 



Dividing the denominator into the numerator and adding one 
to the result one gets 



0.99206 + 0.91^12z"^ + 0.85007z“^ + 



-2 



However j looking at the graph (Fig. 132) the precise 
and the approximated solution start diverging after a while 
This is due to the fact that in 

E = 



s + C 



if C goes to zero the equation becomes 



E = 1/s 



which is the final value that the approximated solution 
approaches . 
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I 



This fact is very important because it leads to these 
Important conclusions: 

1. If goes to zero with time, first check the Laplace 
transform and the initial equation substituting by zero. 
If by Inspection they do not lead to the same type of 
solution, this method does not apply. 

2. Verify, before doing the division if any of the 
varying coefficients, unexpectedly, goes to infinity; by 
dividing by zero, for Instance. 




Figure 132 

X-scale = 1.0 units/inch — Precise solution 

Y-scale = .5 units/inch + Approximate solution 
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IV. CONCLUSIONS 



This technique, for the successfully area covered in 
this paper (linear differential equations with time varying 
coefficients), showed to have a good precision even for 
relatively large iteration steps. The labor necessary to 
complete a problem by handwork is not very Involved and 
needs only a small algorithm for computer solution with 
C.P.U. times much smaller than existing numerical methods. 

However, with time varying coefficients that go to zero 
with time, attention must be paid to the application of 
this method. As a general rule, this method does not apply. 

In some cases, the final equation for the approximate 
solution, before the division is performed, must be altered 
in order to not have time varying coefficients (referred as 
Cj^ in Section IIC) in the numerator. Caution is necessary 
to verify, before starting the division steps if any varying 
coefficient unexpectedly has a value of infinity; by division 
by zero v;hen substituting a particular value of Cj^. 

The technique tested in this paper does not apply to the 
case of linear differential equations with varying coefficients 
only. Linear differential equations with constant coefficients 
and nonlinear differential equations are, as well, solved by 
this technique. 

These areas, in the opinion of the author, are valuable 
as topics for other theses. 
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